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SECTION  1 


INTRODUCTION 


.--Tt>is  report  shows  how  constitutive  equations  for. 
elastic/strain-hardening  materials  are  formulated  for  computer 
programs  which  calculate  wave  propagation.  The  constitutive 
equations  considered  here  are  for  Pile  Driver  granite  and  were 
developed  by  i.  Sandler  of  Paul  Weidlinger  Consulting  Engineer. 
(Reference  1) 

Prior  to  performing  calculations  of  the  Pile  Driver  event, 
preliminary  calculations  were  made  by  several  groups  to  investigate 
whether  there  is  agreement  among  various  computer  codes.  The  pre¬ 
liminary  calculations  involve  wave  propagation  for  which  the  input, 
geometry,  material  properties,  finite  difference  mesh  size  and 
integration  time  step  were  specified  by  the  Defense  Nuclear  Agency 
(DNA)  project  officer.  Although  the  calculations  represent  dynamic 
loading  /<v a  spl i t-Hopki nson  Bar  testing  device,  no  reference  is 
made  to  any  physical  tests.  The  spl it-Hopkinson  Bar  was  thought 
merely  to  be  a  useful  configuration  within  which  codes  can  be  com¬ 
pared.  Comparison  among  the  results  of  the  various  preliminary 
calculations  reveals  some  agreement  and  some  disagreement.  A  study 
conducted  by  Agbabian  Associates  indicated  that  part  of  the  differ¬ 
ences  among  the  results  of  the  preliminary  calculations  was  due  to 
different  methods  of  coding  the  constitutive  equations.  When  each 
of  the  methods  was  suitably  refined,  however,  all  gave  comparahle 
results. 


This  report  documents  the  methods  of  coding  the  capped 
model  which  were  used  during  the  investigation.  The  exercises  which 
were  performed  and  results  which  were  found  to  be  acceptable  are  also 
discussed. 
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The  cap,  which  is  elliptical  and  tangent  to  the  failure 
envelope  at  their  point  of  intersection,  is  given  by 


f2  -  (J.  -  v)2  +  r2(J2  -  «)  -  0 


(2-4) 


in  which 


V  =  L  +  R  X(L)  X * (L) 


(2-5) 


Q  =  (X(L)  3 2{l  +  R2 [X *  (L)]2| . (2-6) 


and 


X  (L)  = 


|fl['  -  ('  - 
(a  +  c 


+  c 


L  <  B 


L  >  B 


(2-7a) 


X 1  (l) 


_  }^(>  - 

1° 


L  <  B 


I.  >  B 


(2- 7b) 


The  hardening  parameter,  L,  is 


"  J[  9(Jr  V3!)/^ dt 


(2-8) 


where 


=  W 


V3! 


V3!- 


(2-9) 


*2  -  w2  4>r . <2-,oi 
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A 

The  material  coefficients  R  and  W  are 
R  =  3.0 
W  =  JO7  ksi 

The  fracture  surface  and  typical  caps  are  illustrated  in  Figure  2-1. 


Bulk  and  Shear  Moduli 

The  bulk  shear  modulus  B  and  shear  modulus  p  are 
assumed  to  be  constant 

B  =  7000  ksi 

U  “  3000  ksi 


Incremental  Stress/Strain  Relationship 

The  incremental  stress  is  related  to  the  elastic  component 
of  strain  by  the  expression 


d°ij  -  (°  -  HKh  *  4eu) . (2'") 


where 


de. .  -  de? . 
|J  'J 


(2-12) 


and  the  plastic  stra i n- i ncrement 


is  defined  by 


d 


P 

e.  . 

'J 


.(2  13) 


k 


"X 
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All  methods  of  coding  the  associated,  plastic  potential 
flow  rule  begin  with  these  statements.  Thus 

do..  =  x(de.  .  -  A  +  2u(de..  -  A  . (2-14) 

'J  \  kk  3 o^f  ij  \  ij  3o.j/ 

or,  in  terms  of  stress  invariants, 

dJ.  -  3K(dckk  - 3A  107) . (2',5a) 

and 

iJJl  -  — !i—  or: .  (del .  -  l4-| . (2-lSb) 

V  2  yjr  u\  U  Sajj/ 


where  primed  quantities  denote  deviator ic  components  of  stress  and 
strain. 

The  final  stress  state,  represented  by 

/o..)  =  ho..)  +  do . (2-16) 

V  ,J/new  V  'j/0,d  <1 

or 

M  "  I'')  ,H  +  dJ<  . l2',7a) 

'  'new  '  'old 

H  -  (V3!),  ‘‘V3 I . <2-'2b> 

'  'new  '  'old 


must  lie  on  the  curve.* 


*The  present  example  is  most  easily  illustrated  with  a  perfectly  plastic 
material.  Modifications  for  a  work  hardening  material  ore  discussed  at 
the  end  of  the  chapter. 
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or 


f(oij) 


=  0 


(2-1 8a) 


(2-1 8b) 


The  mathematical  statement  of  this  constraint,  used  by  all 
contributors  to  the  present  exercise,  \* 


df 


or 


3f 

3  o . . 
U 


0 


df 


3f  ,  .  3f 

dJ  +  - 

3J1  ' 


dWJ‘ 

z 


o 


(2-1 9a) 


(2-1 9b) 


At  this  point,  two  different  approaches  are  used  to  obtain  the  states 
of  stress  and  plastic  strain.  The  first  approach  to  be  described  is 
referred  to  as  Method  A  and  is  used  by  several  groups  engaged  *n  ground 
shock  calculation.  The  subscript  t  ’iutfsfed  to  indicate  a  trial  or 
temporary  state  which  lies  in  the  forbidden  region  outside  f  *  0  as 
indicated  in  Figure  2-2.  The  increment  in  stress  or  stress  invariant 
from  the  old.  equilibrium  state  to  the  trial  state  is  denoted  by 


(do.j )  ,  (djJ  and 


Equations  2-15a  and  •  Jpb  may  be  rewritten  as  follows: 
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FIGURE  2-2.  METHOD  OF  CORRECTING  FIN/'L  STRESS  STATE  FROM  A  TRIAL  STATE 


Substituting  Equations  2-20a  and  -20b  into  2-'$b 


+<)t 


9KAf. 


+  fv 


MyAo!  .a.'  .f_ 
. ^  -*■ 11 

t  2t/]I«JjT 


-  0  .  .  .  (2-21) 


Equation  2-2 1  may  also  be  expressed  by 


9KAf!  +  uAf„  -  f 


n 


(2-22) 


The  expression  on  the  right-hand  side  of  Equation  2-22  is  considered  to  be 
equal  to  the  difference  between  the  trial  state  and  the  surface  f  =  0. 
Thus 


(Jl)t*  (VJ2)J  1  0  -  fl(dJ,)t  +  . 


A  is  then  computed  from  the  expression 


A  = 


H-  NJ 


9KfT2  +  pf. 


(2-24) 


Tl 


The  plastic  strain  increments  can  be  found  by  substituting  this  value 
of  A  into  Equation  2-13-  The  increments  in  Jj  and  can  be 
fouod  by  substituting  A  into  Equations  2-20a  and  -20b.' 

In  an  alternative  approach,  referred  to  as  Method  B, 
Equation  2-14  is  substituted  into  Equation  2-19a  as  follows. 


f . . (Xde 
•J 


kk 


Af  )6..  +  2yf..(de..  -  Af..)  *  0 

11  ij  ij  ij  ij 


(2-25) 
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and  A  is  computed  from  the  expression 


A 


A  (dc 


1 1 

kk'  Hi 
kk^u"" 


+  2yde. . f . . 

,  ,  ,  .JJ.  JJ. 

+  2  (if .  .  f . . 

'J  'J 


(2-26) 


If  the  stress  point  lines  on  the  capped  portion  of  the  yield 
surface,  Equations  2-24  and^-26  must  be  modified  to  account  for 
motion  of  the  cap.  In  thepresent  model,  this  modification  is  as 
follows. 


Method  A 


A  = 


(V31>t 


9Kf  J  +  Uf  2U 


3f 


_2 

3L 


V4 


+  4  +  f  33 


(2-24a) 


Method  B 


>(dEM)f»  ^  2udEi  |f  i  j 


A,UfU  +  2V‘:  ■■  *  9 


'J 


'll 

3L 


Vfn  +  4  + 


33 


(2-26a) 


The  values  of  A  obtained  by  Equations  2-24  and  2-26  are  the  same  if 
the  derivatives  of  f  with  respect  to  stress  and  stress  invariants 
are  the  same.  The  derivatives  differ  if  the  stresses  used  to  evaluate 
them  differ.  Some  groups  use  the  trial  stress  invariants  ( J ^ ) t ,  (/J£) 
to  evaluate  the  derivatives.  Other  groups  use  the  stresses  atthe 
previous  state  (o  j  j )  ^  to  evaluate  the  derivatives.  If  the 
difference  between  the  trial  and  the  old  states  is  large  enough, 

or  if  f  is  a  sufficiently  nonlinear  function  of  its  arguments,  then 
noticeable  differences  between  the  two  approaches  will  develop. 
Appreciable  differences  uo  not  develop  if  either  of  the  above  condi¬ 
tions  is  absent.  In  some  of  the  examples  given  below,  the  trial  and 
the  old  states  are  kept  close  by  the  device  of  subdividing  the  strain 
increment  into  many  smaller  increments. 
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SOME  PROBLEM  AREAS  IN  CODING  MODELS  OF  THIS  TYPE 

The  main  difficulties  in  coding  this  model  arise  because 
important  parameters  of  the  model  are  implicit  function1-  of  stress 
and  plastic  strain  quantities.  Before  the  new  stresses  and  plastic 
strains  can  be  computed,  functions  which  depend  on  them  must  first  be 
evaluated.-  The  dilemma  is  relatively  easily  resolved  when  the 
function  f  corresponds  to  the  fracture  surface  because  the  value* 

'*o7  3f/3ojj  do  not  usually  change  rapidly  with  changes  in  ojj. 

Almost  any  technique  involving  iteration  or  subdivision  of  the  strain 
increments  leads  to  a  correct  answer  so  long  as  enough  steps  are 
taken  from  (oijJojj  to  (°ij)new*  *n  ^act»  the  slope  of  the 
fracture  surface  is  constant,  an  exact  relationship  between  the  trial 
and  fina'  stress  states  can  be  found  in  closed  form  as  follows.  Assum¬ 
ing  the  fracture  criterion  be  given  by 


f,  "  *'J2  "  aJl  "  C  -  0 


(2-27) 


The  associated,  plastic  potential  flow  rule  leads  to 


(2-28) 


yjj7-c  -°Jn 

“  2 

9Ko  +  G 


(2-29) 


The  strain  hardening  or  cap  portion  of  the  model  is  much  more 
difficult  to  code.  The  rate  at  which  convergence  to  the  right  answer 
proceeds  depends  greatly  on  the  technique  of  iteration  or  of  sub- 
divi  ling  the  strain  increment. 


Two  basic  approaches  were  used  for  finding  the  correct  cap 
and  stress  point  in  th<  examples  reported  below.  One  approach  is  to 
rely  completely  on  an  iterative  technique.  Such  an  approach  hopes  to 


accommodate  do, 


which  are  of  the  same  order  as 


hi) 


old 


and 
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correspondingly  large  changes  in  the  cap  parameter.  It  is  argued 
that  this  procedure  is  consistent  with  the  finite  changes  in  stresses 
and  velocities  which  occur  ir.  time-marching  computations  of  shock 
wave  p  ropaga t i on . 

An  alternative  approach  is  to  subdivide  or  split  the 

strain  increment  into  a  series  of  subincrements  in  which  do.,  is 

much  smaller  than  (o:.)  ...  This  has  the  effect  of  enabling-* 

•  j  old 

equations  such  as  Equation  2-19  to  be  regarded  as  linear  difference 
equations.  The  derivatives  3f/3cj.,  etc.  are  still  implicit 
functions  of  plastic  strain,  and  hel^ce  an  iterative  scheme  such  as 
a  Newton-Raphson  technique  or  a  modified  Euler  approach  which 
averages  the  matrix  of  generalized  stress/strain  coefficients  is 
still  required.  It  is  argued  that  this  is  consistent  wit',  the 
assumption  that  the  stress/strain  rules  being  coded  are  continuous 
in  slope  and  function  value.  One  difficulty  wh/efrT joth  methods 
must  overcome  is  illustrated  in  Figure  2-3,  which  shows  how  a  trial 
stress  increment  may  penetrate  both  the  cap  and  the  fract  re  surface. 
The  questions  which  the  code  must  resolve  automatically  are 

a.  Which  surface  ihould  the  stress  point  be  on? 

b.  Which  point  of  the  surface  selected  in  a  should  the 
stress  point  be  on? 

These  questions  are  discussed  further  in  Appendix  A. 
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^APPLICATION  OF  THE  CAPPED  MODEL  IN  TWO  EXERCISES 


The  results  of  two  exercises  are  reported  below.  The  first 
exercise  represents  a  spl i t-Hopkinson  bar  apparatus  with  a  granite 
specimen,  whose  constitutive  relations  are  governed  by  the  capped 
model,  inserted  between  steel  bars.  Figure  3"1  shows  the  dimensions, 
zone  size,  time  step  and  properties.  The  stress  path  followed  by 
elements  in  the  center  of  the  granite  is  illustrated  in  Figure^3-2 . 
Since  most  of  the  zones  undergoing  plasticity  are  in  a  state  of 
approximately  uniaxial  stress  only  a  small  part  of  the  capped  model 
is  exercised  in  the  computation. 

To  expand  the  scope  of  the  present  investigation,  further 
exercises  were  performed  in  which  five  different  stress  paths  are 
followed.  The  paths  are  specified  by  the  combinations  of  principal 
strains  and  are  listed  in  Table  3-1.  The  paths  in  P^ane 

which  an  elastic  material  would  follow  when  subjected 
to  these  strains  are  shown  in  Figure  3“3* 

SPLIT-HOPKINSON  BAR  CALCULATION 

Some  of  the  differences  in  the  spl i t-Hopkinson  bar 
calculation  which  prompted  the  present  study  are  shown  in  Figure  3_i» 
through  3"9-  Group  1  correctly  coded  the  model  using  Method  A  as 
described  above.  Group  2  correctly  coded  the  model  using  Method  B. 
Besides  using  a  different  method  to  compute  A,  the  main  differences 
in  coding  the  material  model  were: 

a.  In  evaluating  functions  of  stresses,  Group  1  used  the 

trial  state  stresses  (ot)  whereas  Group  2  used  the 
previous  state  stresses  (°0 Id)* 

b.  Group  1  subdivided  the  strain  increments  into  smaller 
subincrements  whereas  Group  2  used  the  full  increment. 

Thus,  although  both  methods  are  in  principle  correct,  significant 
differences  occur.  No  exact  solution  is  available,  so  that  it  is 
impossible  to  say  with  assurrance.  that  one  computation  is  better 
than  the  other.  However,  due  to  splitting  the  strain- increment , 
the  results  obtained  by  Group  1  have  a  better  chance  of  being  more 
accurate.  Due  to  performing  a  greater  number  of  operations  in  the 
mater<a^  property  subroutines  of  the  computer  program,  the 
calculation  performed  by  Group  1  is  also  more  expert ve  than  it 
would  be  if  splitting  were  omitted. 
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FIGURE  3~1.  CAP  MODEL  TEST  PROBLEM 


-Jj,  KSI 


ROXIMATELY  UNIAXIAL  STRAIN  UNTIL 
EFACTIONS  FROM  EDGES  ARRIVE 


STREfS  PATH  FOR  AN  ELEMENT  IN  CENTER  OF  GRANITE 
ON  CENTERLINE  IN  SPLIT-HOPKINSON  BAR  COMPUTATION 
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TABLE  3-1.  PRESCRIBED  STRAINS  FOR  STRESS/STRAIN  EXERCISES 


3 

0.179) 


e11 

e22  “  e33  ”  ne11 

2.38  x  10-Z* 

2.38  x  10 

2.38  x  10'3 

2.38  x  10‘3 

2.38  x  10-2 

2.38  x  10‘2 

-4 

-4 

3.2  x  10 

1.13  x  10 

3.2  x  10‘3 

.13  x  10‘3 

-2 

.  -2 

3.2  x  10 

1.13  X  10 

6.4  x  10‘2 

2.26  x  10'2 

-4 

3.5  x  10 

0.626  x  lO-2* 

3.5  x  10'3 

0.626  x  JO'3 

2.8  x  10‘2 

0.501  x  10'2 

5.6  x  10“2 

1.00  x  10'2 

4.13  x  10-4 

*0.875  x  10_1* 

6.61  x  10"3 

-1.4  x  10'3 

1.15  x  lo”2 

-2.45  x  10'3 

3.31  x  10~2 

-7-01  x  10“3 

-4.96  x  10~5 

-3.43  x  10~5 

-4 

-4 

-7.95  x  10 

-5-57  x  JO 

Positive  sign  indicates  compressive  strain 
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(a)  APPROXIMATE  STRESS  PATHS  FOR  STRESS/StfrAi^EXERCISES  (NEW  ORIGIN) 
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S\ 

\  w^'X  _ —  INITIAL  POSITION  OF  THE  CAP 


(b)  FRACTURE  CRITERION,  CAP  AND  ELASTIC  PATHS  (OVERALL  VIEW) 
FIGURE  3-3.  STRESS  PATHS  IN  STRESS/STRAIN  EXERCISES 
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STRESS/STRAIN  EXERCISES 

To  provide  a  basis  for  comparing  various  methods  of  coding 
the  capped  model,  a  series  of  example  problems  involving  prescribed 
strain  paths  has  been  performed.  The  approximate  paths  in  stress 
plane  are  shown  in  Figure  3~3b  and  the  prescribed  strains  are  shown 
in  Table  3*1.  The  stress  states,  based  or  //^dependent  computations 
by  five  groups,  are  ^j&ulated  in  Table  3*2.  Acceptable  variation  is 
considered  to  be  ±\%.  The  amount  of  inelasticity  increases  as  the 
path  number  increases  from  1  to  4.  Thus,  in  path  1  (hydrostatic 
compression)  it  makes  negligible  difference  whether  large  or  small 
strain  increments  are  used,  since  the  response  is  almost  linear. 
However,  in  path  4  (similar  to  unconfined  compression)  it  makes 
considerable  difference  to  some  methods  whether  large  or  small  strain 
increments  are  used,  since  the  response  is  markedly  nonlinear. 


SUBDIVIDING  STRAIN  INCREMENTS  (SPLITTING) 

Due  to  the  extra  cost  required  to  perform  splitting,  which 
can  be  a  substantial  increase  for  finite  difference  codes,  some 
studies  were  performed  to  indicate  the  rate  of  consequence  with 
number  of  splits  or  strain  subdivisions.  This  information  is 
summarized  in  Tables  3*3  through  3*5-  For  the  present  model  and 
magnitude  of  strain  increments  it  appears  that  less  t>sn  10  splits 
of  the  strain  increment  lead  to  acceptable  conveyance.  However, 
even  10  splits  may  impose  an  unacceptable  burden  on  the  cost  of  a 
large  calculation  if  all  strain  increments  for  3ll  zones  are 
split  indiscriminantly.  Thus  it  it  desirable  to  have  an  automatic 
procedure  for  determining  when  splitting  is  necessary  and  how  many 
splits  wi 1 i  lead  to  acceptable  convergence.  At  present,  no 
universally  accepted  criterion  has  been  developed  which  applies  to 
all  situations.  Criteria  have  been  devised  for  individual  problems, 
however.  One  of  these  is  illustrated  in  Figure  3*10. 


ITERATIVE  TECHNIQUES 

Even  if  the  stress  and  strain  increments  are  kept  small, 
by  the  device  of  splitting,  it  is  still  likely  that  a  stress  point 
which  should  lie  on  the  cap  at  thf»  end  of  a  subincrement  does  not 
lie  on  it.  This  happens  because  motion  of  the  cap  depends  on 
plastic  strain  increments  and  is  therefore  not  known  until  the  end 
of  an  increment  or  a  split.  Since  the  motion  of  the  cap  can  be  of 
the  same  order  as  the  stress  increment,  deviation  between  stress 
point  and  cap  may  develop.  To  avoid  o'  minimize  this,  iteration 
is  performed  within  a  strain  increment  or  split.  Two  types  of 
iteration  are  described  in  Appendix  A  which  have  been  used  for  the 
capped  model . 
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TABLE  3-2.  STRESS  STATES  CORRESPONDING  TO 
STRAIN  PATHS  IN  TABLE  3*1 


Path  1 


(E  j  +  e2  +  Ej)/3  (o1  +  02  +  °3^/r3 

(in. /in.)  (ksi) 


0.000233 

5.00 

0.00238 

49-9 

0.0238 

4  98 

Path  2 


E 


1 


( in. /in. ) 


(ksi) 


(in. /in.) 


2  3 

(ksi) 


4.65 

46.4 

462 

900 

0.000113 

0.00113 

0.0113 

0.0226 

3.41 

34.0 

339 

668 

Path  3 

E1 

F  =  E, 

a  „  =  0, 

1 

2  3 

2  3 

(in. /in.) 

(ksi) 

(in. /in.) 

(ksi) 

4.48 

0.0000626 

2.75 

44.7 

0.000626 

27.5 

357 

0.000501 

2.9 

670 

0.01 

425 

Path  4 

E1 

0, 

E„  =  E, 

0 „  -  0- 

1 

2  3 

2  3 

(in. /in. ) 

(ksi) 

(in. /in.) 

(ksi) 

3.67 

-0.0000875 

0.665 

57-5 

-0.0014 

14.0 

99.7 

>  -0.00245 

27.0 

275 

-0.00701 

95.0 

Path  5 

E, 

°1 

E„  =  E, 

0-  =  G« 

1 

2  3 

2  3 

(in. /in.) 

(ksi ) 

(in. /in.) 

(ksi) 

-0.00001*96 

-0.893 

-0.0000348 

-0.805 

-0.000795 

-3.0 

-0.000557 

-3.0 
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TABLE  3-3.  VALUE  OF  HARDENING  PARAMETER  FOR  PATH  4  (WITH  CAP) 

OBTAINED  BY  GROUP  2  FOR  DIFFERENT  AMOUNTS  OF  SPLITTING 

Hardening  Parameter  L,  ksi 


Strain 

E1 

Strain 

Increment 

AE 

UU 10 

AE/100 

AE/1000 

0.0G04n 

0 

0 

0 

0 

C. 00661 

19.2 

58.7 

59 

58.9 

0.0115 

66.2 

109 

110 

110 

0.0331 

409 

352 

339 

359 
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TABLE  3-4.  VALUE  OF  AXIAL  STRESS  (a,)  FOR  PATH  4  (WITH  CAP)  OBTAINED 
WITH  DIFFERENT  AMOUNTS  OF  SPLITTING  BY  THREE  GROUPS. 
DIFFERENCE  BETWEEN  GROUPS  1  AND  2  FOR  AE/IOOO  IS  ATTRIBUTED 
TO  DIFFERENT  DEFINITION  OF  STRAIN  VARIABLES. 


Group 


Strain 


Axial  Stress ,  Oy  ks i 


Strain  Increment 


E1 

AE 

AE/10 

AE/100 

AE/100G 

0.000413 

3.67 

3.67 

3-67 

3.65 

0.00661 

58.7 

57.7 

57.7 

57.7 

0.0115 

88.2 

99-9 

100.0 

99-9 

0.0331 

272 

277 

278 

278 

AE/N  (N< 1 00) 

AE/250 

AC/ 1000 

0. 000413 

3.67 

3-67 

3.67 

0.00661 

57.5 

57.5 

57.5 

0.0115 

99. 

1 

99.2 

99-1 

0.0331 

272. 

273 

273 

AE 

AE/2 

0.000^13 

3.67 

3.67 

0.00661 

57-9 

57.9 

0.0115 

100.7 

100.6 

0.0331 

275.0 

275.0 

TABLE  3-5.  VALUE  OF  RADIAL  STRESS  (o2)  FOR  PATH  A  (WITH  CAP)  OBTAINED 
WITH  DIFFERENT  AMOUNTS  OF  SPLITTING  BY  THREE  GROUPS 


Group  Radial  Stress,  a^,  ksi 


Strain 

Strain  Increment 

E1 

AE 

AE/10 

AE/100 

AE/10Q0 

2 

0-000^13 

0.664 

0.664 

0.662 

0.649 

0.00661 

10.7 

14. 1 

14.1 

14.1 

0.0115 

15-6 

27.3 

27-3 

27.3 

0.0331 

57-9 

96.4 

97-0 

97-0 

AE/H  (N<100) 

AE/250 

AE/1000 

1 

0.000413 

0.665 

0.665 

0.665 

0.00661 

1 4.0 

1 4.0 

14.0 

0.0115 

26.9 

26.9 

26.9 

0.0331 

94.2 

94.4 

94.3 

AE 

AE/2 

3 

0.000413 

0.664 

0.664 

0.00661 

13.9 

13.9 

0.0115 

27.0 

27.0 

0.0331 

94.5 

94.5 
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Main  uncertainties  are: 


Computation  of  Q— Individual  stress  components  or  Individual  strain 
components  or  frequencies  thereof,  such  as  Q  -  |  d  e.  ,  |  +  |  de^  |  +  |  d  | 
are  used.  i 


Choice  of  CRIT — magnitude  of  Q  which  signals  that  a  stress  or 
strain  Increme. .t  'arge  enough  to  require  splitting  has  occurred. 

Choice  of  DELTA — number  of  splits  that  are  going  to  be  allowed,  up 
to  100  in  the  present  example. 


AA4855 


FIGURE  3-10.  AUTOMATIC  SPLITTING 
31 
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A 

SECTION  4 

SUMMARY  AND  CONCLUSIONS 


In  the  present  study,  various  methods  of  ceding 
constitutive  equations  are  examined.  It  is  found  that  certain 
details  of  coding  can  affect  the  accuracy  with  wiich  stresses  are 
computed  from  prescribed  strains  and  strain  incr exeats.  A  series 
of  test  cases,  including  computation  of  wave  propagation  in  a 
spi i t-Hopkinson  bar  apparatus  and  prescribed  strain  paths,  is 
described  which  will  help  reveal  coding  errors  in  cases  where 
complicated  constitutive  equations  are  being  coded  for  the  first 
time. 


Specific  findings  of  the  present  study  include  the 

fol lowing: 

a.  Inaccuracy  in  computing  stress  from  strain 
increments  can  be  reduced  to  an  acceptable  level 
by  the  devices  of  iteration  for  the  values  of 
variables  which  are  implicit  functions  of  their 
arguments,  such  as  Q  in  Figure  2-}r  and  by 
subdividing  or  'splitting'  strain  Increments. 

b.  Iteration  and  splitting  may  increase  the  cost 

of  numerical  computations  significantly.  Hence, 
It  is  necessary  to  have  automatic  controls  which 
introduce  iteration  and  splitting  only  when 
necessary.  Widely  applicable  criteria  for  this 
have  not  been  found.  Some  suggestions  for  the 
specific  model  and  stress  paths  considered  in 
the  present  report  are  given  above. 
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APPENDIX  A 

METHODS  OF  CODING  THE  CAPPED  MODEL 


One  iterative  method  invo1ves"’adjusting  the  stresses  in 
different  wavs  according  to  whether  otr  lands  in  Region  2  or  in 
3  or  4  of  Figure  A-1.  If  otr  lands  in  2,  the  ultimate  stress 
state  is  defined  to  be  on  the  fracture  surface.  A  linear 
approximation  to  the  surface  in  the  vicinity  of  otr  is  made,  and 
a  previously  derived  correction  formula  (which  uses  the  plastic 
potential  flow  rule)  is  applied  to  the  stresses.  The  linear 
approximation  is  necessary,  according  to  one  group,  because  it  is 
otherwise  not  possible  to  work  out  the  correction  in  closed  form. 

If  otr  lands  in  Region  3  or  4,  the  ultimate  stress  state 
is  defined  to  be  on  the  cap.  However,  at  this  state,  the  cap 
corresponding  to  the  ultimate  stress  state  has  not  yet  been  defined. 
Thus  it  is  necessary  to  iterate  for  a  stress  state  which  coincides 
with  a  cap..  First,  a  trial  stress  s£»tfe  is  computed  based  on  the 
position  of  the  cap  at  the  previous  step  and  the  plastic  potential 
flow  rule.  When  a  cap  which  supposedly  corresponds  to  the  new 
stress  state  is  computed,  it  is  found  that  the  cap  and  stress  state 
do  not  coincide  because  the  equations  used  to  relate  the  stress 
state  to  the  cap  parameter  are  based  on  the  assumptions  of  small 
increments.  To  reconcile  the  position  of  the  cap  to  the  stress 
state,  a  binary  search  is  used. 

As  illustrated  in  Figure  A-1,  a  temporary  upper  bound 
point  on  the  next  cap  (S+)  is  defined  by  assuming  the  strain 
increment  to  be  entirely  elastic.  A  temporary  lower  bound  point 
on  the  next  cap  (S“)  is  defined  as  shown.  The  midpoint  between  S" 
and  S+  is  adopted  as  a  trial  state  on  the  next  cap.  A  trial  value 
of  L  and  all  partial  derivatives  are  evaluated  at  the  trial  state. 
From  these  quantities,  a  new  cap  is  computed.  If  the  trial  state 
does  not  lie  on  the  cap  within  some  tolerance,  changes  are  made  in 
the  trial  state  by  an  iterative  procedure  such  that  the  stress  point 
eventually  coincides  with  a  cap.  Similar  procedures  are  used  to 
determine  a  stress  point  on  the  fracture  surface,  which  is  an  easier 
matter  since  it  is  stationary. 

For  finite  element  applications,  it  is  convenient  to  cast 
the  instantaneous  incremental  stress/strain  relations  in  the  fol¬ 
lowing  matrix  form 

do  »  C  dc  (A-1) 


A-1 


The  iterative  scheme  shown  in  Figure  A-2  is  similar  to  the  modified 
Euler  formula.  The  procedure  to  obtain  the  stress  increment  during 
a  step  is  described  below. 


Wk-H  *  2  {  ^k  +  ^k+l}lAElk+1 


(A-2) 


where  (C ) ^  and  [C]|<+i  are  the  instantaneous  stress/strain 
relationships  at  the  beginning  and  the  end  of  the  step,  respec¬ 
tively.  Since  (Cj^+j  depends  on  the  stresses  at  the  end  of  the 
time  step,  successive  corrections  to  [CJk+j  and  to  {AoJk+j 
are  made  until  the  following  function  of  stress  increments 
reaches  a  stationary  value  within  a  prescribed  tolerance: 

E  -  lAoJ  +  |Ao2|  +  |  Aa3|  . (A-;) 


FIGURE  A-2.  ITERATION  SCHEME  TO  EVALUATE  STRESSES 


